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Abstract. We present a structural analysis of the graphene/Ru(0001) system 
obtained by surface x-ray diffraction. The data were fit using Fourier-series expanded 
displacement fields from an ideal bulk structure, plus the application of symmetry 
constraints. The shape of the observed superstructure rods proves a reconstruction 
of the substrate, induced by strong bonding of graphene to ruthenium. Both the 
graphene layer and the underlying substrate are corrugated, with peak-to-peak heights 
of (0.82 ± 0.15) A and (0.19 ± 0.02) A for the graphene and topmost Ru-atomic layer, 
respectively. The Ru-corrugation decays slowly over several monolayers into the bulk. 
The system also exhibits chirality, whereby in-plane rotations of up to 2.0°in those 
regions of the superstructure where the graphene is weakly bound are driven by elastic 
energy minimization. 
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1. Introduction 

Nanostructured materials have attracted increasing interest in recent years, due to 
their potential in practical electronic applications. One of these, graphene, has been 
theoretically investigated since the 1940s [1|. The discovery in 2004 that freestanding 
graphene may be prepared [2] led to an explosion of interest in this material due to 
its unique electronic properties and possible practical utilization |3j. Graphene is a 
single sheet of carbon atoms arranged in a honeycomb structure, which was believed to 
be thermodynamically unstable under ambient conditions, due to the Mermin- Wagner 
theorem Nowadays the stability of graphene is explained by postulating small out-of- 
plane corrugations, leading to lower thermal vibrations [5|[6]. On crystalline substrates 
the formation of superstructures and the concomitant corrugation of graphene provides 
template functionality [7]. The influence of the substrate and the formation of the 
superstructure is believed to change the electronic bandstructure and the electronic 
properties |8j|9], due to bond formation and charge-transfer phenomena (7,10 12 . The 



characterization of the graphene-metal interface structure is of crucial importance, 
because measurements of the electronic transport properties require making metallic 



contacts 13 . Surface x-ray diffraction (SXRD) is a powerful investigative tool for this 
system, since the diffraction intensity is perfectly described in a single scattering picture 
and is unaffected by density-of-states effects or electrostatic forces. 

Graphene grown on transition metals forms single-domain superstructures with 
high degrees of structural perfection [I4[jl8] . Early reports on graphene on Ru(0001) 
proposed a superstructure in which (12 x 12) unit cells of graphene sit on (11 x 11) 
unit cells of ruthenium, ("12-on-ll") [T4lfT5], while other studies proposed an 11-on- 



10 superstructure 16 . However, recent SXRD results showed unambiguously that the 



reconstruction is in fact a surprisingly large 25-on-23 superstructure 17 . A comparative 
study between density functional theory (DFT) calculations and scanning tunneling 
microscopy (STM) experiments showed the structure to be composed of regions of 
alternating weak and strong chemical interactions of the graphene with the Ru-substrate 



19,20 



Here, we detail the atomic structure of the graphene/Ru(0001) system, determined 
with sub-Angstrom resolution from SXRD data. In addition to quantifying the 
corrugation, we also show that the best model exhibits the formation of chiral domains, 
resulting in a lower symmetry (p3) compared to graphite (p3ml). This unexpected 
property may have an important impact on e.g., the use of this system as a template 
for molecular chiral recognition. We argue that this breaking of the symmetry is driven 
by energy minimization based on elastic energy considerations. 

2. Experimental 



Sample preparation and the SXRD measurement setup at the Surface Diffraction station 
of the Materials Science beamline, Swiss Light Source, have already been detailed in 17 . 
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It was demonstrated from simple simulations of the 25/23 superstructure rod (SSR) that 
the substrate must also be corrugated, since oscillations with the appropriate periodicity 
of approximately 1.0 out-of-plane substrate reciprocal lattice units (r.l.u., 2tt/c) on the 
SSRs only start to appear if one includes a corrugation of the substrate. Here, we 
present further SXRD data from the same sample, which in addition to the SSRs now 
also includes in-plane data. 

Because of the very large number of atoms involved in the superstructure, it is 
impossible to fit each atomic position individually. Instead, we have parametrized the 
structural model using a small set of physically reasonable parameters. The in-plane 
and out-of-plane deviations of the atomic positions from an ideal flat structure of the 
graphene and of the uppermost layers of the Ru-substrate are described by a 2D Fourier- 
series expansion. We truncate this series after the fourth Fourier component, since 
higher orders could not be resolved in the diffraction data. The displacement field of 
the system is allowed to adopt the lower p3 symmetry, since this is the lowest symmetry 
still compatible with the apparent measured sixfold diffraction symmetry, which only 
arises because of the superposition of the two possible terminations of the hexagonal 



close-packed (hep) substrate 21 . Because the p3 symmetry allows chiral structures, 
we have to sum over the signals from domains of each enantiomer, and assume a 50 % 
distribution. 

Details of the implementation of the Fourier expansion and of symmetry constraints 
are given in the Appendix. Here, we only discuss those aspects needed to understand the 
results. First, it is important to note that because the 25-on-23 structure contains 2x2 
corrugation periods, only the even Fourier components, that is the second and fourth, 
must be considered. This is also demonstrated by the absence of signal at the 22/23, 
24/23, . . . superstructure rods. For each atom within the supercell, the in-plane and 
out-of-plane deviations Ax, Ay, and Az are described by the two Fourier components. 
In total, both the graphene and ruthenium require nine fitting parameters each in order 
to describe their corrugations. 

In addition to the 18 corrugation parameters we introduce a factor, A, which 
describes an exponential decay of the substrate corrugation amplitude with substrate 
depth z 

A(z)=A -exp(-z/X). (1) 

This decay applies to all the three amplitudes used for the description of the substrate 
displacement function. We fix the minimum distance from the substrate to the graphene 



layer, dc-R U , to 2.0 A 19 22 since our model is relatively insensitive to this parameter 
within physically sensible limits (±0.1 A) . The parameter g?r Ui _r U2 is the distance 
between the first and second Ru-atomic layer. Lastly, a global scaling factor S is 
required, resulting in a total of 21 free fitting parameters. 

We begin by defining regions of the supercell, where we consider a flat graphene 
layer lying commensurably 25-on-23 on top of a flat Ru-substrate (see Fig. [I]). The gray- 
shaded region in Fig. []Jb) indicates where the first of the two C-atoms within a 'normal' 
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Figure 1. (color) (a) Graphene (black) on top of ruthenium (red, first layer) and 
green (second layer). Three regions shaded gray, green, and red, are highlighted, and 
explained in the text, (b) A zoom into the lower left corner of the 25-on-23 supercell. 



graphene unit cell sits on top of a Ru-atom of the topmost substrate atomic layer (red 
atoms), and the second atom sits on top of a Ru-atom from the second substrate atomic 
layer (green atoms, the hep position). Henceforth, we refer to this as the (top, hep) 
region. Using the same arguments, the red area is the (hcp,fcc) region, and the green is 
the (fcc,top) [23]. 



Fitting ||24j was performed using GenX 25 , an optimization program using the 
differential evolution algorithm, which helps avoid getting trapped in local minima |26|. 
The errors on the fitted parameters are estimated by an increase in the goodness-of-fit 
of 5 %. 

We fit g?r Ui „r U2 to the CTR-data alone [Fig. |2](a)] as this is sensitive to small 
differences in the interplanar spacing of the topmost two Ru-atomic layers but is largely 
insensitive to the form of the weakly scattering superstructure. The best fit had an 
R-factor of 5.2 %, for g?r Ui _;r U2 = (2.080 ±0.003) A, which should be compared to a bulk 
value of 2.141 A. This equates to a contraction of 2.8 %, in agreement with 27 -29 



3. Results and discussion 



The starting model for the search of all the other parameters was a strained 25-on-23 flat 
graphene layer lying commensurably on a flat ruthenium bulk structure. The best fit for 
the SSR- and in-plane data has an R-factor of 13.4 % [Fig. [2|b) and (c)]. The peak-to- 
peak corrugation height of the graphene is (0.82 ±0.15) A, in agreement with 



14 



16 



30 



whereas that of the uppermost Ru-atomic layer is (0.19 ± 0.02) A and is out of phase 
with respect to the graphene corrugation (Fig. |3|. The exponential decay length of 
A = (7.0 ±0.4) A means there is still approximately a tenth of the distortion of the first 
Ru-atomic layer at a depth of four Ru-atomic layers. This strongly supports the idea of a 




Figure 2. (color online) (a) The (1,0)-CTR. Only the scaling factor and rfR Ul -Ru 2 
were used to fit the data, (b) Fit of the two superstructure rods, (c) In-plane map 
of the superstructure reflections around the (1,0)-CTR position at I = 0.4 r.l.u. The 
areas of the circles are proportional to the scattering intensities. 



chemisorbed graphene layer with significant interaction with the substrate |7[ [22j[3l}|33] . 

Details of the final structure are summarized in Fig. |ij Figure Qa) shows a 
clear corrugation of the graphene with the hills lying in the weakly bound (hcp,fcc)- 
region. The hills have a triangular shape, in remarkable agreement with earlier STM 



data 14 15 . Although in-plane movements of up to (0.25 ±0.03) A of the graphene are 
observed [Fig. |4^b)], the bond lengths are distorted by less than 0.1 A. This requires a 
twisting motion and indeed the in-plane movements exhibit a chiral signature, in which 
the largest movements occur at the steepest flanks of the hills, as one might expect, 
based on simple elastic strain considerations. Note that this feature emerged naturally 
from the fitting and was not implemented a priori into the model. The biggest rotation 
angle of the hexagons is 2.0° found on the flanks as well as on top of the hills. 
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The elastic energy was calculated to test the physical validity of the presented 
parametrization approach and the resulting model. It takes into account the in-plane 
and out-of-plane displacements of surface atoms from their 'ideal' positions due to the 



25/23 surface reconstruction. From our model, we calculate an elastic energy 34-36 



due to strain of 9.3 eV per supercell, assuming zero strain for a flat 25-on-23 graphene 
layer 37 . Fitting the data to the higher p3ml-symmetry results in an increase in the 
elastic energy by 83 %, while the R-factor of 14.7 % is significantly higher than that for 
the p3-symmetry. Even if we were to assume zero strain for a flat graphene layer having 
the bulk graphite in-plane lattice constant [Fig. |4^c)], this has no significant influence 
on the energy difference between the two different symmetry models. A histogram of all 
the bond lengths in the graphene superstructure demonstrates that the implementation 
of the lower p3-symmetry allows the bond lengths to be more preserved relative to bulk 
graphite. 



We are aware of an independent low-energy electron-diffraction (LEED) study 38 



on the same Ru single crystal using the same graphene preparation and characterization, 
where the authors claim a corrugation of the graphene layer of 1.5 A, and a corrugation 
of the topmost ruthenium layer of 0.23 A. In that study the system is described by a 
p3ml-symmetry and the unit cell is cut down to one of the four inequivalent sub-unit 
cells in order to reduce computational time. The SXRD simulation of the coordinates 
extracted from the LEED-study led to a significantly higher R-factor of 34.0 %. The 
reason for the discrepancies, which are outside the error bars, are not yet resolved, 
although possible explanations are the already-mentioned restriction to p3ml-symmetry 
and a 12-on-ll superstructure - a full dynamical scattering LEED calculation of the 
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Figure 4. (color) (a) Top view of the counterclockwise twisting enantiomer resulting 
from the fitting procedure: The graphene shows the lowest lying atoms to be in 
(top,hcp)-region, whereas the hill maxima appear in the (hcp,fcc)-regions. Clear 
triangular-shaped hills are observed, (b) The in-plane displacements of the same 
enantiomer, magnified by a factor of 10, from the ideal bulk positions. The distortions 
are largest on the flanks of the hills, (c) Histogram of the bond lengths in the graphene 
layer. The model with p3 symmetry allows the carbon hexagons to twist, and most 
of the bonds are stretched by less than ±0.04 A compared to the bulk bond length of 
graphite (1.421 A) (blue, dashed line) or a flat 25/23 superstructure bond length 
(aR U /\/3 x 23/25 = 1.4373 A) ( green, dotted line). Enforcing the higher p3ml- 
symmetry causes larger distortions in the bond lengths. 
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system with p3-symmetry is presently beyond computational capabilities. In addition, 
the fact that LEED only probes the topmost layers, while SXRD demonstrates that 
significant vertical displacements occur down to four atomic layers of the Ru-substrate, 
might also play an important role. 

4. Summary and conclusion 

In summary, we have determined the graphene/Ru(0001) structure in unsurpassed 
detail. This was only possible by adopting a parametric Fourier description of the 
superstructure using only a small number of physically reasonable parameters. Up to the 
mirror-symmetry breaking the final model agrees excellently with previous STM studies. 
We find a graphene and ruthenium corrugation peak-to-peak height of (0.82 ± 0.15) A 
and (0.19 ± 0.02) A, respectively. The ruthenium corrugation is out of phase with that 
of the graphene and decays exponentially down to a depth of several ruthenium layers. 
Importantly, we have also discovered the new and potentially highly significant property 
of areal chirality in the in-plane movements, which are most evident on the flanks of 
the hills of the corrugation. We propose that this symmetry-breaking phenomenon is 
induced by elastic energy minimization of the graphene layer. To test the validity of 
this, we calculated the elastic energy of the graphene superstructure to be 9.3 eV, less 
than two thirds of that for the p3ml case. 
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Figure Al. (color) Different corrugation shapes generated by different out-of-plane 
phase values. The blue regions are the strongly bound "valleys", the red highlighted 
regions show the weakly bound "hills" . 



Appendix 

In the following, the implementation of the symmetry constraints and the Fourier 
expansion to the graphene-on-ruthenium model will be briefly described. The 
displacement dr of an atom sitting at point r is expressed by its 2-dimensional Fourier 
series 

dr 1 = £ Kl >t • sin[2vr( S x + ty) + # j (A.l) 



K,t = V<* 2 + Bi/, 4>\,t = arctan^/^) (A.2) 

where s, t e {0,2,4} are the orders, A\ t , B l st are the Fourier coefficients, <p\ t are the 
phases of the corrugation and i G {x,y,z}. Note that the phase of the out-of-plane 
displacements influences the valley and hill shapes and positions of the corrugation 



allowed by the p3 symmetry (Fig. Al). 
Since 

sin(/ + 0) = sin(/) cos(0) + cos(/) sin(0) (A. 3) 

and equating A 1 and B l to 

A i = K l ■ cos(^) 

B i = K l -sin^), (A.4) 
one can rewrite Eq. |A.1| as 

d^ = A U ■ sin[2vr(sx + ty)] + B* >t ■ cos[2vr(sx + ty)]. (A.5) 

s,t 

The rotation operators used for the description of the p3-symmetry are Ri and R2, 
which in a hexagonal coordinate system describe a 120° rotation counterclockwise and 
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dr[R(r)] 
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► 



dr(r) = R-^drpEl^)]} 
= R{dr[R- 1 (r)]} 



Figure A2. The p3-symmetry constraint operators. R is defined as a rotation by 
120° counterclockwise around the origin, while R _1 is the rotation clockwise by 120° 
around the origin 



clockwise around the origin, respectively (Fig,A2), are given by 
.I -l)' R2 = ( -I ) ■ 



Ri 



(A.6) 



It can be easily shown that R4 = R 2 1 = R. dr has to fulfill the p3-symmetry 
constraint, which results in 

R _1 {dr[R(r)]} = R{dr[R _1 (r)]} = dr(r). (A.7) 

The relations in Table Al follow by inserting Eq. A. 5 and Eq. |A.6 in Eq. A.7 



A U = 


A z 




A z 


= 

^i,t = 


AX 

A t-(s+t) 
AX 

A -(s+t),s 


+ A t,-(s+t) 

- A y 


= ~ A -(s+t),s 
AX 

- A t,-(s+t) 


B s,t = 


jdz 

n t,-(s+t) 




jdz 

n -(s+t),s 


JDX 

Bit = 


n t-(s+t) 

JDX 


+ B t-(s+t) 
n -(s+t),s 


- - B v 

TJX 

- D t,-{s+t) 



Table Al. The relations of the Fourier coefficients A\ t and B\ t (i £ {x, y, z}). 
Regarding the considered Fourier components in the analysis, the zeroth order is the 
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23/23 reflection, the first and third order components correspond to the 24/23 and 26/23 
systematic absences, respectively, and the 25/23 to the second order component. Hence 
the fourth order refers to the 27/23- reflection, and along the /i-direction (equivalent to 
the fc-direction) one can limit (s,t) = (2, 0) and (s,t) = (4, 0). For the sake of simplicity, 
we describe here only the implementation of to the second order. 

From Eq. A. 5 and Table Al , one can derive the following expressions for the single 



components of dr which describe the displacement field. We do not include the orders 
(s,t) for the sake of simplicity 

dr z = A z ■ sin(2vr • 2x) + A z ■ sin[27r(-2y)] 
+ A z ■ sin[27r(-2x + 2y)\ 
+ B z ■ cos(2tt • 2x) + B z ■ cos[27r(-2y)] 

+ B z ■ cos[27r(-2z + 2y)\ (A.8) 
dr x = A x ■ sin(27r ■ 2x) - A y ■ sin[27r(-2x + 2y)\ 
+ (A y - A x ) ■ sin[27r • {-2y)\ 
+ B x ■ cos(2tt ■ 2x) - B y ■ cos[27r(-2x + 2y)] 

+ {B y - B x ) ■ cos[2tt • {-2y)\ (A.9) 
dr y = A y ■ sin(27r ■ 2x) - A x ■ sin[27r • (— 2tt ■ 2y)\ 
+ (A x - A y ) ■ sin[2vr(-2x + 2y)\ 
+ B y ■ cos(2tt ■ 2x) - B x ■ cos[2tt • (-2tt ■ 2y)} 

+ (B x - B y ) ■cos[27i(-2x + 2y)}. (A.10) 

So, we obtain six fitting parameters for the displacement field dr, namely A x , 
A y , A z , B x , B y , B z . Since we lock the phase for the fourth order to be the same as 
that for the second order, there will be nine fitting parameters. The effect the fourth 
order harmonic with a locked phase has on the structure is shown in Fig. A3 For an 



amplitude up to 0.25 of that of the second harmonic, the low regions in the structure 
will be flattened out. 
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Figure A3. Effect of the implementation of the fourth harmonic: fi represents the 
second harmonic, / 2 the fourth. Their sum for an amplitude of the fourth harmonic 
up to 0.25 of that of the second harmonic makes the low regions flatter. 
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